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The accuracy of the many-body perturbation theory GW formalism to calculate electron-phonon 
coupling matrix elements has been recently demonstrated in the case of a few important systems. 

However, the related computational costs are high and thus represent strong limitations to its 
widespread application. In the present study, we explore two less demanding alternatives for the 
calculation of electron-phonon coupling matrix elements on the many-body perturbation theory 
level. Namely, we test the accuracy of the static Coulomb-hole plus screened-exchange (COHSEX) 
approximation and further of the constant screening approach, where variations of the screened 
Coulomb potential W upon small changes of the atomic positions along the vibrational eigenmodes 
are neglected. We find this latter approximation to be the most reliable, whereas the static COHSEX 
ansatz leads to substantial errors. Our conclusions are validated in a few paradigmatic cases: 
diamond, graphene and the Ceo fullerene. These findings open the way for combining the present 
many-body perturbation approach with efficient linear-response theories. 

PACS numbers: 71.15.Qe,71.38.-k 


I. INTRODUCTION 

Electron-phonon coupling (EPC) occupies a prominent 
role in various fields of condensed matter physics, includ¬ 
ing phonon-mediated superconductivity, photoemission 
band gap renormalization, current carriers inelastic scat¬ 
tering or the lifetime of hot electrons. Concerning the cal¬ 
culation of EPC matrix elements on the ab initio level, up 
to now, mainly density functional theory (DFT) 1 and its 
perturbative linear response extensions (DFPT(P® 

have 

been applied, providing important information at the mi¬ 
croscopic level. 

Recently, several studies questioned the accuracy of the 
DFT-based EPC matrix elements when obtained with 
(semi)local functionals such as LDA or PBE.® By way 
of example, the electron-phonon coupling involv ing spe¬ 
cific electronic and phonon modes in graphene 5 6 and 
graphite,® the value of the electron-phonon coupling po¬ 
tential to states at the Femi level in the electron-doped 
fullerenep® superconducting bismuthates and transition- 
metal chloronitrides,® or the renormalization of the pho¬ 
toemission band structure of pentacene 11 and diamond 
crystals,® were shown to be affected by a significant un¬ 
derestimation of the EPC matrix elements when calcu¬ 
lated within DFT and (semi)local functionals. 

As a cure to such problems, many-body perturba¬ 
tion theory (MBPT) techniques within the so-called GW 
approximation 13 *® showed a clear impr ovem ent when 
compared to available experimental data.® 2 ®- 12 Unfor¬ 
tunately, the kind of theories that are available within 
DFT, and in particular the powerful DFPT formalism,®^ 


are not yet available within the framework of MBPT. 
Existing GW calculations of EPC matrix elements have 
been therefore based on the frozen-phonon approach, ne¬ 
cessitating a stepwise displacement of the atoms along 
the phonon modes with an explicit evaluation of the elec¬ 
tronic structure on the GW level for each step. This de¬ 
manding approach cannot be reconciled with very large 
unit cells and is thus only feasible for zone-center or 
zone-boundary phonon modes. Going beyond the frozen- 
phonon approach, in order to access EPC matrix ele¬ 
ments on the GW level for a large number of electron 
and phonon wave vectors, remains a considerable chal¬ 
lenge to the ab initio community. Significant work is still 
needed to come up with a scheme that may allow to gen¬ 
eralize the pioneering above-mentioned studies to a much 
larger set of systems. 

In the present study, we explore the merits 
of the st atic co ulomb-hole plus screened exchange 
(COHSEX)l 14 l 15 l 1 M 2 tl approximation to the GW self¬ 
energy for the calculation of the electron-phonon cou¬ 
pling matrix elements in the isolated fullerene molecule, 
in diamond and in graphene. Our definitions of the rele¬ 
vant EPC matrix elements follow previous studies rela ted 
to the superconducting transition in the fullerides^® 
the zero-point motion renormalization of the gap in 
diamond 12 -* and the gap opening through electron-phonon 
coupling in graphene,® respectively. For the same set of 
systems, we also explore the accuracy of the constant¬ 
screening approximation, namely the assumption that 
the screened Coulomb potential can be considered, to lin¬ 
ear order, as a constant upon small changes of the atomic 
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positions. Whereas COHSEX leads to non-negligible dis¬ 
crepancies, this latter approximation will be shown to be 
robust and accurate. The present results offer impor¬ 
tant perspectives in making the GW formalism amenable 
to the study of condensed matter phenomena involving 
electron-phonon coupling. 


II. METHODOLOGY 

We briefly introduce the many-body perturbation the¬ 
ory Green’s function formalism, providing a solid frame¬ 
work for the calculation of quasiparticle energies E. In 
such an approach, the one-body quasiparticle eigenvalue 
equation reads: 

(zZ! + ^-( r /) ^(r) + V H (r)0(r) 

+ J dr'£(r,r';£)0(rO=i^(r), 

where V lon and V H are the ionic and Hartree poten¬ 
tial, respectively. The self-energy £(r, r';E) replaces 
the well-known exchange-correlation potential of den¬ 
sity functional theory or the exchange operator in the 
Hartree-Fock formalism. In general, it is non-local, 
energy-dependent and non-Hermitian. Derived within 
Schwinger’s functional derivative approach to perturba¬ 
tion theory,® the GW approximation to the self-energy 
leads to: 


Yp w ( y , y * \E) = J due i5u} W(r,r'-uj)G{Y,Y'-E + uj), 

where G and W are the time-ordered one-particle Green’s 
function and the dynamically screened Coulomb poten¬ 
tial, respectively, and S = 0 + a small positive infinitesi¬ 
mal. 

The Coulomb-hole plus screened-exchange (COHSEX) 
representation of the GW self-energy was originally intro¬ 
duced using a time representation of G, W and £) 14 l 2 ° l 21 l 
We follow Hybertsen and Louic 15: by using the following 
spectral representations of G and W: 


G(y,y';E + co) 
W (r, r 7 ; ui) 


= £ 
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E + uj — e n — z<5sgn(£ n — fi) 
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where G has been written in terms of one-body eigen¬ 
states <j> n and eigenenergies £„, typically starting DFT 
Kohn-Sham or Hartree-Fock solutions, and W in terms 
of its spectral function B. These expressions allow to 
obtain the pole structure of both G and W. From the 
residue theorem, one then rapidly obtains an exact de¬ 
composition of £: 


occ 

Y, SEX (r, r'; E) = - ^ W(y,y’\E — £„)^(r)^*(r') 

n 

E coh (y,y';E) = )^n(r')^ [ du:' 

^7 Jo E — £ n — UJ 

where V indicates the principal value. T, SEX , which in¬ 
volves a summation over the occupied states only, orig¬ 
inates from the poles of G. It is called the screened ex¬ 
change interaction in analogy to the bare exchange term 
that can be obtained by replacing W with the energy- 
independent bare Coulomb potential v. £, COH originates 
from the poles of W and represents the Coulomb-hole 
contribution, since it can be shown to be related to the 
interaction of an electron with its related adiabatically 
built correlation hole. 

The static approximation to the exact COHSEX de¬ 
composition assumes that (E—e n ) ~ 0 for all (n), leading 
to simplified static screened exchange and Coulomb-hole 
expressions, namely: 

OCC 

£f t f£(r,r';0) = ^l¥(r,r';0)i(r)C(r') 

n 

ESS&r.r'iO) = ^£<Ar t ( r )^( r, )^( r , r, ;°) 

n 

= *<Kr-r , )W(r,r';0). 

Here, W = {W—v) is the difference between the screened 
and bare Coulomb potential. As a result, besides being 
a static approximation, this Coulomb-hole term is also 
local in space. Such a static COHSEX approximation 
(labeled COHSEX here below) was shown to yield too 
large gaps in the case of semiconductors.® By way of 
example, in the present case of the Ceo molecule, the 
COHSEX gap is found to be 5.3 eV, i.e. 0.4 eV larger 
than the ~4.9 eV experimental gap.®- 4 ' Nevertheless, this 
has to be compared to the starting 1.6 eV DFT-LDA 
Kohn-Sham gap which is dramatically too small. 

While it cannot be claimed that the static COHSEX 
approach is a good approximation to absolute quasipar¬ 
ticle energies, we emphasize that we are interested in 
quasiparticle energy differences upon small (infinitesi¬ 
mal) atomic lattice motions. The main assumption on 
which we rely to calculate the electron-phonon coupling 
within the COHSEX approximation is that the varia¬ 
tions of the dynamical contribution to the self-energy 
can be neglected. This can be a priori rationalized by 
emphasizing that dynamical interactions are driven by 
the plasmons dynamics, collective excitations less sensi¬ 
tive to small atomic displacements than single-particle 
excitation energies and wave functions. It remains, how¬ 
ever, that besides the static approximation, the spatially 
local character of the static COH term is at odds with 
the nonlocality of the full GW self-energy. 
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In order to further justify the following results concern¬ 
ing the constant screening approach, namely the second 
approximation we explore in this study, it is instructive 
to consider the GW plus Bethe-Salpeter formalism!^®] 
The latter is a many-body perturbation theory approach 
concerned with describing the linear response of a sys¬ 
tem with respect to a time-dependent external perturba¬ 
tion and thus the Bethe-Salpeter equation (BSE) is the 
MBPT analogue to time-dependent DFT. At the heart 
of the GIT/BSE approach lies the variation (dGW/d\), 
where the “perturbation” A is the one-body Green’s func¬ 
tion G. The most common approximation that has been 
shown to be remarkably accurate 31 " 33 ' is to replace the 
GW self-energy by its static COHSEX approximation, 
and to consider further that (dW/d A) = 0. It is such a 
simplified scheme we aim to explore in the present study, 
differing in the fact that the perturbation A is now in¬ 
duced by a vibrational distortion of the system. 


III. TECHNICAL DETAILS 

The many-body GW and COHSEX calculations for 
the isolated Ceo fullerene are performed using the Fiesta 
code, an implementation of the GW formalism within 
a Gaussian basisP^^H We start from DFT-LDA eigen¬ 
states calculated with the Siesta package 37 and a triple- 
zeta plus polarization (TZP) basis®! for the description 
of the valence orbitals combined with standard norm- 
conserving pseudopotentials. 391 We use the resolution of 
the identity technique (RI-SVS) with an even-tempered 
auxiliary basis set formed of four Gaussians per angular 
channel (up to d character) with exponents between 0.2 
and 3.2 Bohr -2 . With such a basis, our GoUo@LDA 
ionization potential and HOMO-LUMO gap are found 
to be 7.33 eV and 4.39 eV (B3LYP geometry), respec¬ 
tively, starting from a 1.66 eV DFT-LDA gap. This is in 
close agreement with the 7.45 eV and 4.40 eV plane wave 
GoRo@LDA values of Ref. 0H1 (LDA geometry). The 
GoHo@LDA gap remains smaller than the ~4.9 eV ex¬ 
perimental value, an issue that will be addressed by the 
partially self-consistent scheme used in the present study, 
as discussed below. 

Following Ref. [8] we use the relaxed structure and 
phonon eigenmodes generated within a DFT-B3LYP ap¬ 
proach using a 6-311G(d) basis. This was shown to pro¬ 
vide vibrational frequ encie s in excellent agreement with 
Raman measurementsConcerning the calculation of 
the EPC matrix elements, our present study follows 
the approach described in Refs. (8191 where the electron- 
phonon coupling strength associated with the (electron- 
doped) lowest unoccupied 3-fold degenerate ti u molecu¬ 
lar orbital (LUMO) has been explored within DFT and 
GW. The specific choice of the LUMO level is dictated 
by the physics of the phono n-m ediated superconducting 
transition in the fullerides ! 42 ! 43 ! We adopt the standard 
definition of the effective electron-phonon coupling po¬ 
tential V ep that enters e.g. the McMillan formula for 


the superconducting transition temperature through the 
dimensionless parameter A = N(Ep)V ep , where N{Ep) 
is the density of states at the Fermi level. Namely, V ep 
reads in the so-called molecular limit jMEU 


v* =—y — y 

9 M ^ ^ 

v u m—1 

Here, M is the mass of the carbon ions, is the fre¬ 
quency of the vibrational mode with index (v) and u" is 
the vibrational polarization vector. As already discussed 
before, the electronic states with energy (e m ) are lim¬ 
ited to the 3-fold degenerate LUMO level. From group 
theory analysis it follows that only the H g and A g vi¬ 
brational modes can couple to these states. A schematic 
representation of the ti u level splitting as a function of 
the deformation amplitude along a H g mode is provided 
in Fig. [Tfa). The EPC matrix elements are consequently 
related to the slopes (<9e m /<9u''), showing a strong depen¬ 
dence on the formalism adopted to calculate the energies 
e m . These energy derivatives are calculated within the 
frozen phonon approach, using a symmetric five points 
finite-difference formula. 

Such an effective electron-phonon potential was calcu¬ 
lated for Ggo hr Refs. 15151 both within DFT using diverse 
(semi)local and hybrid exchange-correlation functionals 
and GW. While DFT-LDA values were shown to signif¬ 
icantly underestimate the coupling energy, DFT calcula¬ 
tions performed with global hybrids containing 20%-30% 
of exact exchange and the GW approach were found ^ to 
compare favorably to gas phase experimental dataP§I 47 ! 
As a drawback of DFT calculations with hybrids, it was 
evidenced that the resulting V ep potential would quickly 
increase with the amount of exact exchange, yielding the 
usual question of the proper choice of the functional pa¬ 
rameters for a given system. 

In the present study, we are concerned with explor¬ 
ing simplified GW schemes, i.e approximations reducing 
the needed computational effort. As such, quasiparti¬ 
cle energies are evaluated at a single-shot COHSEX and 
a partially self-consistent evCOHSEX level, where the 
quasiparticle energies are reinjected self-consistently in 
the construction of G and W, while keeping the start¬ 
ing Kohn-Sham LDA wave functions unchanged. Our 
reference point are partially self-consistent e\GW calcu¬ 
lations with self-consistency on the eigenvalues as pre¬ 
sented in Ref. [9] for the calculation of the EPC matrix 
elements in Cq o- Such a simplified approach to full self- 
consistency, justified by the dramatically too small start¬ 
ing point DFT-LDA gap for small isolated molecules, 
has been shown in the case of gas phase organic systems 
to produce mo re accur ate ionization potentials and elec¬ 
tronic affinities,! 34 ! 35 ! 48 ! ^! together with imp rove d optical 
excitation energies at the GW /BSE level ja 6 A 2 H54] w j ien 
LDA or PBE Kohn-Sham eigenstates are chosen as start¬ 
ing points. Consistently, when starting from DFT-LDA, 
EPC matrix elements in Cq q were also found to be in 


de n 


d \T 
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FIG. 1: Schematic representation of changes in the electronic structure under phonon distortion: a) in Cqo, the 3-fold LUMO 
levels (blue) split under a distortion (red) along a H g phonon mode, b) In graphene, an energy gap A E at the Dirac point is 
opened (red) when distorting the equilibrium structure (blue) along the r — A] or K — E? g phonon modes, c) In diamond, the 
lowest lying conduction bands at F /X (blue) are renormalised by their coupling to a zone-boundary phonon mode. For clarity, 
only the states showing the largest shift are plotted (red dot, see text). 


closer agreement with experiment within evGW as com¬ 
pared to GqWqW 

For the periodic diamond and graphene systems, we 
use the Yambo code which implements many-body per¬ 
turbation theory within a plane wave formalism. 5 ^ The 
starting DFT Kohn-Sham eigenstates are generated with 
the Abinit codeP^ In the case of graphene, the running 
parameters are identical to those used in a previous GW 
study of the electron-phonon coupling in this material^ 
namely we use a = 2.46 A as lattice parameter, 15 a.u. 
distance between the layersp 36 x 36 x 1 k-points for 
the GW calculations in the Brillouin zone (BZ) corre¬ 
sponding to the primitive (2-atoms) cell and the nearest 
equivalent k-grid in the BZ corresponding to the super¬ 
cell needed to describe zone-boundary (q = K) phonons 
in the frozen-phonon approach. The energy cutoff on the 
plane wave basis was set to 60 Ry. All conduction bands 
within an energy range of 45 eV above the Fermi level 
were used to build G and W, and a cutoff of 4 Ha for the 
construction of the dielectric function. The Godby-Needs 
plasmon pole model was employed for the dynamical di¬ 
electric constant entering in W. 

As schematically described in Fig. [ljb), we focus on 
the splitting A E of the degenerate occupied and unoc¬ 
cupied levels at the Dirac point caused by the coupling 
to the qv = T-F^g optical mode at zone center together 
with the zone-boundary qv = K-Aj^ phonon, for which 
the strongest renormalization of the coupling constant 
have been observed upon replacing the DFT-LDA ap¬ 
proach by the Go Wo formalism PS The EPC matrix ele¬ 
ments of interest are thus here related to the derivatives: 
(dAE /du q „), namely (see Ref. [5j): 


(Dl„) 


1 (dA E\ 2 
Aq„ [du^qj 


( 2 ) 


expressed in {eV/A) 2 , where A q „ equals 16 and 8 for 
qv = T-E 2 g and qv = K-A 1 , respectively. We per¬ 
formed “single-shot” Go Wo and COHSEX calculations on 
top of DFT-LDA results, namely the standard perturba¬ 


tive (non self-consistent) GW approach to obtain quasi¬ 
particle energies. 

For diamond, following Antonius and coworkersp- 2 -! we 
have used a relaxed DFT-LDA geometry (lattice parame¬ 
ter a = 3.591 A) with a plane wave cutoff energy of 80 Ry 
and norm-conserving pseudopotentials. All calculations 
were performed with a 8x4x4 uniform Monkhorst-Pack 
grid to describe the supercell necessary for the frozen- 
phonon simulation of the (q = X) phonon. We have 
used 160 bands and a cutoff of 8 Ha for the construc¬ 
tion of the dielectric function and the GW self-energy. 
Following Refs. fT2l and f58l where the diamond band-gap 
renormalization by the zero-point motion was studied, 
we define the relevant EPC matrix elements as the sec¬ 
ond derivative of the electronic eigenvalues with respect 
to the phonon displacements at equilibrium: 


demk _ Ti d 2 

dn vq 2Mu) U q 3u2 q £mk ’ 

where n vq is the occupation of the phonon state (v) of 
wave vector q. The EPC matrix elements are now re¬ 
lated to the temperature dependent renormalization of 
the electronic bands. We evaluate the second-derivative 
of the eigenvalues using a symmetric five-point finite- 
difference formula with displacements of 0.01, 0.02 and 
0.03 A. In Ref. [52] it was already shown that stan¬ 
dard DFT and DFPT methods with (semi)local function¬ 
als strongly underestimate (by ca. 30%) the zero-point 
renormalization of the direct gap of diamond compared 
to experiment, while Ref. [12] demonstrated that Go Wo 
calculations restored the good agreement between theory 
and experiment. We focus here below on the renormaliza¬ 
tion of the lowest conduction band edge at (mk = Xi c ) 
and (mk = ri 5c ) by the coupling to the (vq = X 4 ) 
phonon modes at the zone-boundary (see Fig. [lj:). It 
is for such a phonon wave vector that the EPC matrix 
elements were found to be large and that also the GW 
correction was shown to be the most significant P^O Fol¬ 
lowing Ref. [12] the provided EPC matrix elements are 
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LDA | GpWp |GoWo(W)|COHSEX|COHSEX(WJ 

Graphene 


r -E 2g 

44.3 

68.1 (65.2) 

72.0 

81.5 

83.2 

K-A', 

88.6 

201(187) 

207 

260 

256 


Diamond 


Ti5 c 

1.281 

0.611 

0.657 

1.311 

1.445 

x lc 

-2.009 

-1.031 

-1.066 

-2.124 

-2.273 


TABLE I: Calculated EPC-related quantities (see text) for graphene and diamond. In the case of graphene, these correspond 
to the splitting of the degenerate highest valence and lowest conduction band at the Dirac point under the influence of a 
(q v = T-E 2 g) and (qis = K-A x ) phonon mode and are expressed in (eV / A ) 2 . Number in parenthesis are obtained with a 
plasmon pole energy of 7 eV, close to the graphene 7r-plasmon energy, instead of the 27 eV originally used in Ref. [5] (see text). 
For diamond, these are the results of equation |3| for the Xi c and Tisc conduction states and for the X4 zone-boundary phonon 
mode only (values in eV). This corresponds to twice the zero-point renormalization arising from this given phonon mode (see 
Ref. |12|). 


averaged over the degenerate electronic manifold consid¬ 
ered (namely the 3 degenerate levels at zone-center and 
the 2 degenerate levels at the X lc zone boundary). 

IV. RESULTS 
A. Reference GW calculations 

We first reproduce for validation and for reference the 
previously published GW calculations related to EPC 
matrix elements in graphene, diamond and Cq 0 - Our re¬ 
sults are compiled in Table |T] for diamond and graphene 
and in Table [IT] for Cq 0 , respectively. Our LDA and Go Wo 
values for graphene are close to those found in Ref. 0 
namely e.g. 201 eV/A (G 0 Wo value, present study) for 
the largest matrix element with the (K-A-J phonon, to 
be compared to 193 eV/A in the previous study.The 
difference can be explained by an increased five points, 
instead of only two, finite-difference formula, and by the 
increase of the dielectric matrix size energy cutoff, from 2 
to 4 Ha cutoff. We also verified that changing the Godby- 
Needs plasmon model input finite frequency, from the 
(default) 27 eV value in Ref. [5] to a 7 eV value closer to 
the 7r-plasmon resonance in graphene, does not change 
significantly the calculated EPC strength (see numbers 
in parenthesis in the Table [I]). Such differences are neg¬ 
ligible with respect to the more than 100% increase as 
compared to the LDA value. 

For diamond, the EPC contribution associated with 
the mk = Ti5 C conduction states are in good agree¬ 
ment with the LDA and Go Wo results of Antonius and 
coworkers pJ i.e. we find all meV difference with re¬ 
spect to their EPC in the case of LDA, while we find a 
difference of 8 meV for the Go Wo results.^ These small 
variations can be ascribed to the different pseudopoten¬ 
tials (Ref. [59] showed that pseudopotentials can lead to 
errors up to 50 meV), and to the different convergence 
parameters. 

Since Ref. [T2l focused on the renormalization of the di¬ 
rect band gap at the zone-center, comparison with the 
present results for the phonon coupling to the Xi c elec¬ 


tronic state cannot be made. Our findings are, however, 
very consistent with what is observed for the Tisc state, 
namely a dramatic decrease of the coupling strength with 
the X 4 phonon mode which, again, is the dominant cou¬ 
pling mode. The coupling energy is indeed found to be 
reduced by 52% and 49% from LDA to GW for the ri 5c 
and X lc states, respectively. Since the X lc state is closer 
to the true conduction band minimum, this is a strong 
indication that zero-point motion renormalization of the 
indirect gap will certainly also be strongly affected by 
the GW correction. The full study of such an effect is 
beyond the scope of the present paper. 

For G 6 o, the total evGW coupling potential is within 
8% of that found by Ref. [5] as a result of the larger TZP 
basis used herepSl Our evGW V ep potential is found to 
increase by 44% with respect to the corresponding LDA 
value, to be compared with the 48% increase obtained 
in Ref. [5] with a smaller basis. This indicates the good 
convergence of the GW correction to the LDA value. We 
note that while DFT and GW calculations are known 
to converge very differently with respect to basis size, we 
compare in the present study approaches which are much 
more similar, namely “standard and approximated” GW 
calculations, indicating certainly an even better conver¬ 
gence of the differences we are interested in. 


B. The COHSEX approximation 

We now explore approximations to the full GW cal¬ 
culations performed in Refs. [5F9I12I and reproduced here 
above. We start by the case of diamond. As compared to 
the GoWo calculations, the static COHSEX approxima¬ 
tion is shown to induce errors larger than 100% for the 
relevant electron-phonon coupling energy to the T 15c and 
X lc states. In both cases, the static COHSEX approach 
strongly overestimates the electron-phonon coupling. A 
representation of the errors as compared to the GoW) 
calculations is provided in Fig. [2] 

The case of graphene is of specific interest since, in 
great contrast to diamond and Cqq, the fundamental gap 
reduces to zero in the equilibrium geometry and opens 
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Mode 


LDA 


Theory Experiment 

evGW | COHSEX | evCOHSEX | evCOHSEX(lI / ) ~ 

Fullerene 


H g ( 1) 

4.7 

5.85 

4.99 

6.3 

6.3 

Hg{ 2) 

9.3 

10.1 

9.99 

10.3 

11.9 

Hg{ 3) 

8.8 

12.95 

11.96 

13.7 

15.3 

Hg( 4) 

4.2 

5.5 

4.9 

5.4 

5.7 

Hg{ 5) 

3.98 

4.9 

4.5 

5.3 

5.5 

Hg{ 6) 

1.8 

1.99 

2.05 

2.2 

2.4 

Hg{ 7) 

15.8 

25.4 

23.3 

27.9 

26.5 

Hg( 8) 

13.1 

18.3 

17.1 

19.8 

19.5 

A g (l) 

1.2 

1.9 

1.85 

1.8 

1.2 

Ag( 2) 

7.2 

13.98 

13.2 

15.4 

12.2 

Total A g 

8.4 

15.9 

15.0 

17.2 

13.4 

Total H g 

61.7 

84.96 

78.8 

90.9 

93.1 

Total 

70.2 

100.9 (108.6“) 

93.8 

108.1 

106.6 


TABLE II: Mode-resolved V„ v coupling potential associated with the 3-fold degenerate LUMO and the corresponding V ep 
totals (in meV). For group symmetry reasons, only the H g and A g vibrational modes couple to this electronic state. 
a Ref. 0 

b Ref. Eg] Table V, 
c Ref. \M 
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FIG. 2: Relative error of the different approximations with 
respect to the GoWq results for the “thermal renormaliza¬ 
tion” EPC energies (Eqn. [3| for the lowest conduction band 
of diamond. 


FIG. 3: Relative error of the different approximations with 
respect to the GqWq results for the electron-phonon coupling 
in graphene. 


upon phonon distortions. For the optical phonon mode 
at the zone-center, the COHSEX approximation leads 
again to an increase as compared to the GqWq value, 
even though smaller than in the diamond case with a 20% 
error. Such an error becomes larger for the K-A 1 mode 
with a 29% increase of the EPC-related gap-opening rate. 
In the present graphene case, the error induced by the 
COHSEX approach remains smaller than the 34% and 
56% reduction observed upon using DFT-LDA instead 
of Go Wo- However, the close to 30% error observed for 
the COHSEX calculation of (D^) is clearly significant 
(see Fig. [3]), questioning again the applicability of the 
static COHSEX formalism for EPC matrix elements cal¬ 
culations. 

We now consider the Goo case, where we compare 
evGW and evCOHSEX calculations. A global increase 
of about 7% from evGW to its static approximation is 
observed. Besides the H g ( 4) and A g {l) modes showing 
very small couplings, the static COHSEX approximation 


is observed again to systematically increase the coupling 
strength as compared to evGWp^ This result is consis¬ 
tent with what was observed for diamond and graphene, 
even though clearly the error appears to vary from one 
system to another. 


C. The constant screening approximation 

Following the above-mentioned analogy with the 
Bethe-Salpeter formalism, we finally test another ap¬ 
proximation, that is the constant screening approach. 
Namely, we explicitly calculate the screened Coulomb 
potential for the undistorted structure and upon chang¬ 
ing the positions of the ions along the (R q „) vibrational 
eigenmodes, we keep W frozen to the undistorted value. 
Formally, this amounts to assuming that: dGW/i 9u qt/ ss 
(9G/9u qi/ )W. Such an approach will be labeled GW{W) 
or COHSEX(W) in the following. 
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FIG. 4: Relative error of the different approximations with 
respect to the evGW results for the Tig-related V ep of equa¬ 
tion 0- For convenience, a comparison to the experimental 
data is provided. 

In the case of the fullerene, where calculations are per¬ 
formed within a Gaussian basis, special care must be 
taken in the implementation of such a constant-screening 
approximation. In the auxiliary basis, or resolution of the 
identity approach, the bare and screened Coulomb poten¬ 
tials are expressed in terms of an atom-centered auxiliary 
basis with the following relations: 

[W]p,p‘ = J J drdr'l3{Y)W{r,r')P'(r'), (4) 

W(r, r') = £ /3(r) (S-'mS- 1 )^, /?'(r'), (5) 

AA 

where S is the overlap matrix in the auxiliary basis. Us¬ 
ing now the notation: W and (3 for the screened Coulomb 
potential and the auxiliary basis for the slightly dis¬ 
torted system, the assumption: W(r,r') ~ W(r,r') 
leads straightforwardly to the condition: 

[EQ/3,/3' — Spp\W]p i piSp>pi, (6) 

where Spp =< (3\j3 > is an overlap matrix between the 
auxiliary bases for the perturbed and unperturbed sys¬ 
tems, respectively. 

For Cqq again, we test the constant screening approx¬ 
imation at the COHSEX level only, namely comparing 
evCOHSEX and evCOHSEX(W) calculations of the V ep 
energies. While the plasmon-pole approach used for dia¬ 
mond and graphene is based on a fixed finite frequency 
pole value independent of the geometry, the real axis 
poles contribution to the correlation energy in the con¬ 
tour deformation approach (see Refs. 1191.5111 implemented 
in the Fiesta package changes from one structure to an¬ 
other, leading to difficulties when trying to implement the 
constant-W approach within the frozen-phonon scheme. 

The results of the constant-screening approximation 
are compiled in Tables [I] and [TT] For Ceo, an excel¬ 
lent agreement is obtained within evCOHSEX(W) com¬ 
pared to the corresponding evCOHSEX calculations. 


Comparing the total coupling, evCOHSEX(W) and ev¬ 
COHSEX agree within 1.5 %. 

We now explore the constant-screening approximation 
for periodic diamond. Results are again compiled in Ta¬ 
ble [I] For diamond, we observe a 3.4% and 7.5% dif¬ 
ference between G 0 W 0 (W) and G 0 W 0 for the X lc and 
ri 5c states, respectively. At the COHSEX level, the 
constant screening approximation yields errors of the or¬ 
der of 7.0% and 10.2%, respectively. In summary, while 
the COHSEX approximation was failing for diamond, 
the constant-W approximation leads again to very rea- 
sonnable results, with errors well within 10% at the GW 
level. As illustrated in Fig. [2j this error is much smaller 
than the error induced by the LDA approximation. 

We finally address the case of graphene. Concern¬ 
ing the coupling with the zone-center optical mode, the 
constant-W approach leads to a 5.7% and 2.1% error 
respectively when comparing GoWo(W) to Go Wo and 
COHSEX(W) to COHSEX. In the case of coupling with 
the zone-corner K-A x phonon mode, the discrepancy is 
of 3% and 1.5% applying the constant-screening approx¬ 
imation to GW and COHSEX approaches, respectively. 
Notice that the results for the constant-screening of the 
K-A x phonon are affected by a numerical error of about 
4%, therefore the discrepancy could be slightly larger. 
Anyway this difference is much smaller than the error in¬ 
duced by the standard DFT-LDA calculations, but also 
much smaller than the error induced by the static COH¬ 
SEX formalism as compared to the reference GW calcu¬ 
lation. 


V. DISCUSSION 

Clearly, among the two approximations explored here 
above, the constant-screening stands as a much better 
approach than the static COHSEX for the calculation 
of electron-phonon coupling matrix elements within the 
GW formalism. Overall, the largest error induced by 
the static COHSEX approximation is larger than 100% 
in the diamond case, while it is reduced to 7% in the 
case of the constant-screening approximation to full GW 
calculations. 

The static COHSEX approximation is known to al¬ 
ways overestimate band gaps. As proposed in Ref. El a 
too large band gap should lead to an underscreening of 
the electron-phonon interaction and consequently to en¬ 
hanced EPC matrix elements. This is consistent with the 
fact that, within Hartree-Fock or with increasing exact 
exchange in hybrid functionals, the EPC matrix elements 
are found to steeply increase. The Hartree-Fock (Up) 
and (D|) coupling constants in graphene were found to 
be about 5 and 30 times larger, respectively, as com¬ 
pared to their GW analog,^ while in Geo increasing the 
percentage of exact exchange from 20% to 30% in hy¬ 
brid functionals was found to enhance the V ep energy 
by about 15%P Such an interpretation matches the ob¬ 
servation that in graphene and Geo, the COHSEX EPC 


























coupling constants are larger than their corresponding 
GW reference. In the case of diamond, where we are 
concerned with the second-order derivatives, the effect of 
assuming a static approximation may be more difficult 
to interpret. 

An important observation further is that the Coulomb- 
hole term is a local potential in the static COHSEX ap¬ 
proximation, washing out spatial local-held effects, a ob¬ 
servation that may potentially allow to understand the 
failure of the COHSEX approximation to reproduce ac¬ 
curately the evolution of the band structure with phonon 
deformation. In any case, the large variations of the in¬ 
duced errors from diamond to graphene and Cq o is still 
to be understood, and we will just stand here on the ob¬ 
servation that the static COHSEX approximation cannot 
be trusted to improve on DFT-LDA calculations. 

Concerning the neglect of the gradient of W with 
respect to the ionic positions, such an approximation 
was already tested by Ismail-Beigi and Louie in the 
context of excited state ionic forces within the Bethe- 
Salpeter formalism,^ showing good agreement with ex¬ 
plicit finite-difference ’’exact” BSE calculations for small 
CO and NH 3 molecules. It is commonly assumed that 
the GW /BSE approach is much more resistant to ap¬ 
proximations on W as compared to the GW approach 
for (charged) excitations. This is due to cancellations 
of errors between the electron-electron and electron-hole 
interactions. Namely, any error introduced in W is ex¬ 
pected to affect excitonic interactions and quasiparticle 
gaps in opposite ways. Clearly, the present GW study of 
the variations of a given quasiparticle energy with respect 
to ionic positions cannot benefit from such cancellation of 
errors. Still, the constant-screening approximation turns 
to be a reliable approach to save on computational cost. 

An important consequence of the present findings is 
that once the screened Coulomb potential W( r, r';w) is 
built for the equilibrium geometry, the calculation of 
the variations of the quasiparticle energies with respect 
to the perturbation (A) only requires the evaluation of 
the variations of the Green’s function G with respect to 
the perturbation. This can be performed within stan¬ 
dard DFPT techniques, at least in the case of non-self- 
consistent Go Wo calculations where the Green’s function 
assumes an explicit form in function of the input DFT 
eigenstates. This may invite, for systems such as Cqq, to 
use Go Wo calculations starting from DFT eigenstates ob¬ 
tained with hybrid functionals, which have been shown to 
be a better starting point for finite-size systems 


VI. CONCLUSIONS 


We have explored two approximations for calculating 
self-energy gradients of interest for electron-phonon cou¬ 
pling, namely the two simplifications commonly used in 
the GW/Bethe-Salpeter calculations, that is the static 
COHSEX and the constant screening approximations. 
We explored these approaches in the case of diamond, 
graphene and G 60 . Our findings suggest that the COH¬ 
SEX approximation cannot be trusted to improve on the 
DFT-LDA values as clearly illustrated in particular in 
the case of diamond. On the contrary, the constant 
screening hypothesis, namely assuming that W remains 
to first order constant with respect to small ionic dis¬ 
placements, seems to be reliable, with a discrepancy no 
larger than 10 %, even in the difficult case of graphene 
where the phonon perturbation dramatically affects the 
Dirac cone and the (semi)metallic nature of the graphene 
sheet. Even though a deeper understanding of the ori¬ 
gin of the specific difficulties uncountered by the static 
COHSEX approximation would allow to better rational¬ 
ize the validity and limits of the tested approximations, 
the present results offer promising perspectives to carry 
on such many-body evaluations of the electron-phonon 
coupling gradients with much reduced computer cost on 
realistic systems, including the study of periodic systems 
with arbitrary wave vector perturbation. Further studies 
are however required on a larger set of systems and physi¬ 
cal observables in order to better assess the interest, with 
respect to common DFT calculations, of using the GW 
approach, and its various approximations, for calculating 
self-energy gradients with respect to ionic positions. 
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